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1.  Introduction 


The  behavior  of  an  exploding  metal  initiator  (a.k.a.  exploding  foil  initiator,  or  slapper)  has  been 
studied  extensively  for  single  fuze  detonator  applications,  but  future  combat  systems  (FCS)  will 
require  warheads  with  multi-point  initiators  that  can  produce  simultaneous  detonations.  One  can 
imagine  the  assembly  of  several  resistor/conductor  bridges  of  identical  dimensions  in  an  attempt 
to  achieve  the  desired  goal  of  simultaneous  initiations,  thus  creating  a  maximum  possible 
explosive  force.  However,  variations  of  1%  in  the  resistivity,  thickness  uniformity, 
dimensions/volume,  surface  roughness  and  grain  structures  in  the  bridges  can  cause  more  than 
50%  variation  in  the  initiator  functionality,  resulting  in  a  very  weak  explosive  force.  The 
researchers  at  ARL  have  discovered  a  unique  MEMS  technology  to  prototype  multiple  Cu 
bridges  (1 ,  2)  with  less  than  0.1%  variation  in  the  above  mentioned  parameters,  thus  creating  the 
first  multipoint  fuzes  that  initiated  with  a  high  degree  of  simultaneity.  The  cost  of  these  fuzes, 
when  produced  in  quantity,  is  approximately  one  hundred  times  less  than  that  of  conventional 
fuzes. 

This  technology  has  several  distinct  advantages  over  conventional  discrete  technology.  Typical 
discrete  processes  use  thick-film  screening  printing,  which  limits  bridge  dimension  accuracy  to 
100  pm  and  alignment  resolution  to  25  pm,  and  in  turn  reduces  the  desired  explosive  force.  The 
use  of  semiconductor  fabrication  methods  and  MEMS-enabling  process  technologies  yields  an 
approximately  1000  times  improvement  in  bridge  dimension  accuracy  (0.1  pm  tolerance)  and 
repeatability  of  those  dimensions  across  any  number  of  points.  This  precision  helps  insure  a  high 
current-density  flow  path  across  a  geometrically  accurate  conductor  that  will  be  shared  equally 
amongst  the  points,  contributing  to  a  very  high  degree  of  simultaneity  across  an  array. 

The  fuze  community  benefits  from  many  years  of  experience  modeling  exploding  foil  initiators 
in  terms  of  lumped  parameter  circuits,  with  parameters  based  primarily  on  empirical  data  fits. 
Several  researchers  have  also  employed  finite  element  methods  to  predict  the  spatial  temperature 
distribution  in  a  metal  bridge  wire  resulting  from  Joule  heating,  with  model  results  generally 
valid  up  to  the  metal  melting  point  ( 4  through  7).  Pre-melt  studies  are  important  because 
simulations  predict  the  metal  to  be  in  its  solid  phase  for  90%  of  the  time  before  reaching  the 
experimentally-observed  burst  and  departure  of  a  flyer  or  slapper  (5,  6).  Thus,  it  seems  that  the 
pre-melt  temperature  profile  largely  detennines  the  post-melt  and  burst  behavior,  and  the 
experimentally  measured  kinetic  action  of  the  flyer  material. 

To  understand  how  our  design  or  fabrication  processes  changes  affect  bridge  wire  perfonnance 
and  ultimately  simultaneity,  we  have  simulated  the  behavior  of  an  exploding  metal  bridge.  In 
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this  report  we  present  results  following  the  methodology  of  others  by  focusing  on  the  thermal- 
electric  problem  of  the  pre-melt  phase.  We  have  modeled  the  problem  using  material  properties 
from  two  metals,  Cu  and  Au.  In  future  reports  we  plan  to  address  post-melting  phenomena  as 
discussed  in  section  4. 


2.  Methods 


2,1  Theory 

The  thermal-electric  coupled  problem  is  given  by  the  equation  from  transient  heat  conduction, 

<v.^=°A|2+v.(*-vr)  (i) 

ot 

and  the  equation  for  the  electric  field  E , 

^  =  E-Vcre+aeVE  =  0,  (2) 

ot 

where  the  electrical  conductivity  <je ,  thennal  conductivity  k  ,  and  specific  heat  cY  are  all 
functions  of  temperature.  The  quantity  8pe  /  dt  relates  to  the  diffusion  of  electrons  in  the  metal, 
which  is  very  short  (5)  and  can  be  set  to  zero  in  equation  2. 

It  is  often  advantageous  to  non-dimensionalize  equations  prior  to  applying  finite  difference 
methods,  in  order  to  1)  avoid  rounding  errors  in  taking  finite  differences  of  very  small  numbers, 
and  2)  to  aid  in  presenting  results  and  applying  them  to  a  more  general  set  of  experimental 
conditions.  Non-dimensionalization  involves  choosing  quantities  for  the  fundamental  units  in 
the  governing  equations,  which  in  equations  1  and  2  are  time  t ,  length  L  ,  mass  m  ,  temperature 
T ,  and  charge  C .  Other  quantities  such  as  current  /  ,  thennal  conductivity  k  ,  electrical 
conductivity  ae ,  specific  heat  cv  (at  constant  volume),  electric  field  E ,  or  voltage  V  are  related 

through  these  units  as  follows: 


I  =  C/t 

(3) 

k  =  mL  /  t'T 

(4) 

<je  =  C2t  /  mE 

(5) 

cv=L2lt2T 

(6) 

E  =  mL  /  t2C 

(V) 

2 


V  =  mL2  /  rC  (8 

In  order  to  generate  non-dimensional  variables  for  equations  land  2  with  values  on  the  order  of 
unity,  we  define: 

t'=t/tT  (9 


L*  ~L1  a 


T*  =[T-Tx\/[Tmelt-Tx\ 


C  =  C/I0t, 


PP 

m  =  m  /  <  p — 


where  is  the  period  of  a  sine  wave  for  which  the  first  half  cycle  represents  the  electrical  input 
current  pulse  of  amplitude  70,  a  is  the  bridge  wire  width,  [pm  ( f  )]  is  the  mass  density  of  the 
bridge  wire  metal  at  ambient  temperature  Tx ,  and  fmelt  is  the  melting  temperature  of  the  metal. 
Using  these  definitions,  we  deduce  the  following  relationships: 

_*  f  d  7  d  ~  d  A  [  1 1  ,,  . 

dx  dy  dz  J  \_  a  _ 


Pm  =  Pm  ! 


k  =kt 


Cv  CJU[Tmelt-Tx\ 


E*  =  E  / 


[>.&>)} 


3 


(20) 


V  =v  / 


1 y 

1oLi 


[o-.(r„)]« 

I*0  =  sin 


(21) 


Substitution  of  equations  14  through  19  into  equations  1  and  2  yields  non-dimensional  versions 
of  the  governing  equations, 


CVPm 


8T 

dt* 


=  a\  E*  +V*.(ifcW) 


(22) 


and 


E  »V  cr,  -Ter, V  is  =  0  . 


(23) 


2.2  Finite  Element  Implementation 

Figure  1  shows  the  two-dimensional  geometry  implemented  in  ANSYS  release  1 1  with  a  half 
bridge  wire  width  of  0.5,  and  a  length  of  1 .  The  width  at  the  inlet  and  outlet  regions  of  the  “bow- 
tie”  was  5,  the  angle  of  the  inlet  and  outlet  converging  or  diverging  boundaries  was  45°,  and  the 
fillet  radius  applied  at  the  corners  to  the  entrance  of  the  bridge  was  0.008.  We  modeled  half  the 
actual  geometry  to  conserve  computing  resources,  and  applied  the  boundary  conditions  shown  in 
figure  1. 
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Figure  1.  Finite  element  mesh,  and  boundary  conditions  of  the  fuze  initiator  geometry 
modeled.  The  symbol  V  refers  to  the  gradient  in  the  direction  normal 
to  the  boundary. 

We  modeled  three  geometry  and  material  combinations  as  shown  in  table  1.  Cases  1  and  2 
corresponded  to  those  in  (5,  6),  and  case  3  corresponded  to  a  single  bridge  foil  of  a  multipoint 
wafer  microfabricated  at  ARL.  The  current  pulse  amplitude  for  case  3  was  3800  A,  which  we 
assumed  to  be  equally  distributed  among  15  bridge  wires  of  the  multi-point  experiment  shown  in 
figure  2.  Thus,  the  input  current  amplitude  to  a  single  wire  was  253.3  A. 


Table  1.  Material  and  geometry  combinations  modeled. 


Case 

Material 

Metal  Thickness 
(lim) 

Bridge  Length 
(gm) 

Input  current  pulse 

Amplitude  (A) 

Period  (ns) 

1 

Cu 

1 

1270 

1200 

588.2 

2 

Au 

1 

1270 

1200 

588.2 

3 

Cu 

2.5 

100 

253.3 

370 

5 


Figure  2.  Experimental  current  pulse  and  streak  photo  from  a  multipoint  initiation  test 
performed  at  ARL. 

We  obtained  thermal  material  properties  for  Cu  and  Au  from  (5).  Figures  3  and  4  show  the 
functional  relationships  for  thennal  conductivity  k  and  specific  heat  cv  based  on  temperature. 

We  obtained  electrical  conductivity  in  units  of  1  /  Qm  as  a  function  of  temperature  T  in  units  of 
Kelvin  from  (3),  so  that 

<Je  |Cu  =  |— 4. 12  x  10-10  +  0. 1 13  x  10~5  (T I  \  1604)' 145  J  '  (24) 

for  Cu,  and 

<j e  |Au  =  j— 4.95  x  10-10  +  0. 170  x  10-5  (T I  \  1604)1 178|  '  (25) 


for  Au. 

One  important  modification  to  the  model  inputs  accounted  for  the  fact  that  ANSYS  assumes  any 
2-D  thermal  electrical  elements  to  have  a  thickness  dl  of  unity.  The  current  density  resulting 

from  a  current  input  I0  to  a  line  of  2-D  nodes  having  a  length  a  is  I0  /  adx .  The  actual  3-D  case 
with  current  I2  has  finite  depth  d2 ,  so  that  the  actual  current  density  is  /,  /  ad2 .  Equating 
current  densities  for  the  two  cases  yields  70  =  I2(dl  / a) (a / d2) .  We  readily  recognize  the 
quantity  (<fj  / a)  as  the  non-dimensionalized  depth  d* ,  which  is  always  unity  in  the  ANSYS 
program.  Thus,  the  current  input  /0  to  the  ANSYS  program  for  the  real  world,  3-D  input  I2 
must  be  scaled  by  the  (o  /  d2 )  ratio: 


6 


a 


(26) 


/  =  — -I 
10  ,  2  ' 

d2 


These  arguments  are  not  necessary  for  a  3-D  model,  but  are  for  a  2-D  implementation  as  we  have 
done  for  this  report. 


Figure  3.  Thermal  conductivity  (k)  and  specific  heat  (cp) 
of  Cu  as  a  function  of  absolute  temperature. 

Solid  lines  and  equations  in  the  form  of  y  =  f  ( X ) 
represent  linear  fits  based  on  the  least  squares  method. 


Figure  4.  Thermal  conductivity  (k)  and  specific  heat  (cp)  of 
Au  as  a  function  of  absolute  temperature.  Solid 
lines  and  equations  in  the  form  of  y  =  f  (x) 
represent  linear  fits  based  on  the  least  squares 
method. 
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The  transient  solver  of  ANSYS  attempted  to  automatically  detennine  a  time  step  size,  but  to 
control  solution  stability  and  to  implement  the  current  input  waveform  accurately  we  imposed  a 
certain  number  of  steps  at  which  to  force  a  solution  (called  “load  steps”  in  the  ANSYS  program). 
First  we  estimated  the  time  to  melt  by  considering  a  lumped  parameter  approximation  with 
constant  material  properties.  With  all  spatial  derivatives  set  equal  to  zero,  equations  22  and  23 
reduce  to  a  single  equation 


(27) 


Noting  that  E*  =  /*  /  a  d*<j*e ,  we  set  E*  =  /*  given  our  non-dimensional  scheme.  Thus,  a  quasi¬ 
static  solution  to  equation  27  is 


T*  = 


(28) 


Equation  2 1  gives  the  relationship  for  I*0  ,  and  for  small  t  we  approximated  equation  2 1  as 


Substitution  of  equation  29  into  28  gives, 


T*  « 


(27rti\: 

*  *  * 

CvPm 


(29) 


(30) 


Thus,  the  estimated  time  for  T 


to  reach  the  a  melting  temperature  of  1  is 


* 


f  *  *  \l/3 

c  p 

vr  m 

J 


(31) 
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Figure  5.  Non-dimensional  temperature  of  the  middle  node  versus  non-dimensional 
time  for  the  three  cases  in  table  1.  The  melting  temperature  corresponds 
to  T  =  1 ,  and  time  is  scaled  by  the  input  current  pulse  period  values  in  table  1. 

We  divided  the  quantity  in  equation  31  by  30  for  the  load  step  time  in  ANSYS,  to  insure 
approximately  30  load  steps  at  which  the  program  would  calculate  a  solution  leading  up  to  the 
melting  point  of  the  center  node  of  the  bridge  wire.  Once  a  load  step  was  reached  where  T*  >  1 
for  this  node,  we  halted  the  solution  and  used  a  linear  interpolation  with  the  previous  load  step  to 
approximate  the  exact  time  at  which  the  temperature  of  this  node  exceeded  the  melting 
temperature. 
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3.  Model  Results 


3.1  Transient  Solutions 

Figure  5  shows  the  temperature  distribution  versus  time  for  the  three  cases  in  table  1 .  The  shapes 
of  the  curves  differ  as  melting  time  increased,  with  an  increasing  change  in  curvature  near  the 
end.  A  reason  for  this  trend  is  the  fact  that,  in  the  first  quarter  period  of  the  input  current  pulse 
(0  <  t  <  0.25  ),  the  current  and  therefore  rate  of  heating  monotonically  increased.  During  the 
second  quarter  period  of  the  input  sine  wave  ( 0.25  <  t  <  0.5  ),  the  current  magnitude  dropped, 
and  the  rate  of  heating  logically  followed  the  same  trend. 

Figure  6  shows  the  spatial  temperature  distribution  at  times  that  correspond  approximately  to 
case  3  central  node  temperatures  of  0.3,  0.6,  and  0.9.  As  in  similar  modeling  studies,  the  current 
density  increased  greatly  near  the  comers  as  a  result  of  high  electric  fields,  and  therefore,  the 
highest  rate  of  heating  occurred  in  the  corners.  The  temperature  then  quickly  spread  throughout 
the  rest  of  the  bridge  region. 

The  melting  time  of  case  1  was  within  1.5%  of  the  predicted  time  in  (5)  for  a  similar  set  of 
geometry  and  material.  This  result  increased  our  confidence  that  the  model  was  implemented 
correctly,  or  at  least  as  correctly  as  another  previously  published  study.  However,  the  melting 
time  of  case  2  was  39%  longer  than  the  time  predicted  in  (5),  and  the  melting  time  of  case  3  was 
22%  longer  than  the  experimentally-observed  burst  point.  The  actual  melting  time  would  be 
very  difficult  to  detennine  experimentally  for  case  3  or  any  other  cases  due  to  the  small  sizes  and 
short  time  scales  involved. 
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Figure  6.  Temperature  distribution  of  a  Cu  bridge  foil  for  case  3  at  a)  t*=0.24,  b)  t*=0.31,  and  c)  t*=0.36.  The 
model  results  are  reflected  about  the  line  of  symmetry  shown  in  figure  1 . 

To  probe  further  into  why  the  predicted  value  of  Au  differed  so  much  from  (5),  figure  7  shows 
the  current  density  near  one  of  the  corners.  The  vector  pattern  differed  greatly  from  a  plot  at  an 
equivalent  time  in  (5),  such  that  in  figure  7  all  vectors  appear  to  have  a  negative  or  zero  vertical 
component.  In  (5)  some  discussion  is  given  about  conductivity  gradients  near  the  corners  giving 
rise  to  a  constriction  of  the  current  stream  right  before  melting,  effectively  increasing  the  corner 
radius  of  curvature.  Current  density  vector  plots  show  significant  narrowing  and  then  widening 
of  the  current  stream,  with  significant  positive  current  density  vector  components  in  the  vertical 
direction.  The  effect  is  most  pronounced  for  Au,  but  is  present  for  all  metals  investigated  in  (5), 
including  Cu. 
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Figure  7.  Current  density  at  t*=0.20  for  case  2,  right  before  the  center  of  sample  reached 
the  melting  temperature. 


These  effects  are  not  present  in  (6),  and  thus,  although  not  explicitly  stated  in  (5),  the  effects 
likely  arise  from  a  different  electrical  conduction  model  for  temperatures  above  the  melting 
point.  A  significant  drop  in  conductivity  beyond  the  melting  temperature  would  produce  large 
conductivity  gradients  and  thus  the  vector  plots  shown  in  (5).  In  our  model,  we  only  extrapolate 
the  conduction  relationship  given  in  equation  24  or  25,  which  when  plotted  is  a  nearly  linear 
relationship  changing  by  less  than  a  factor  of  six  over  the  temperature  range  0  <  T*  <  1 .  Such 
gradual  changes  lead  to  relatively  small  conductivity  gradients  for  the  temperature  distributions 
shown  in  figure  6,  whereas  gradients  would  be  much  larger  with  the  empirically-observed  order- 
of-magnitude  changes  in  conductivity  between  solid  and  liquid  phases. 

3.2  Sensitivity  Study 

In  order  to  investigate  the  effects  of  changing  geometry  and  material  parameters  on  the  pre-melt 
behavior,  we  conducted  a  series  of  simulations  for  case  3  where  we  changed  the  bridge  thickness 
d2 ,  bridge  length  L,  comer  radius  R,  taper  angle  6 ,  thermal  conductivity  k,  mass  density  pm , 

specific  heat  cv,  electrical  conductivity  <je ,  and  melting  temperature  Tmelt  each  by  5%.  For  those 

properties  that  varied  with  temperature,  we  applied  the  5%  change  to  the  rate  of  change  factor 
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(i.e.,  the  linear  slopes  in  figure  3  for  cv  and  k).  We  compared  the  time  for  the  central  node  to 
reach  the  melting  temperature  to  the  original  time  to  melt,  and  figure  8  shows  the  results. 
Interestingly,  geometric  factors  such  as  taper  angle  and  comer  radius  did  not  significantly  affect 
the  results,  which  we  expected  based  on  an  early  study  of  potential  distributions  (P).  The 
thickness  was  expected  to  be  a  major  factor  because  the  thickness  directly  affects  the  applied 
current  in  equation  26.  All  material  properties  affected  the  time  to  melt  except  k,  which  was 
surprising  until  we  considered  the  arguments  below. 


Figure  8.  Percent  change  in  time  to  melt  resulting  from  5%  changes  in  individual  geometric 
and  material  property  variables. 


The  one-dimensional  form  of  equation  22  with  constant  material  properties  is 


cvpm  dT*  1 12  ,  d2T 


r  sc 


+- 


ox 


*2 


(32) 


A  solution  to  the  homogeneous  solution  (neglecting  the  forcing  term  I*2 )  by  separation  of 
variables  yields 


f  k* 

t'\ 

v'T  (  *  / 

.  \  .  /  *  / 

exp 

*  * 

n  2 

C.  COS  \X  /  , 

l  +  C,  sin  (x  /  / 

l) 

V  CvPm 

A  ) 

1  \ 

/  z  \ 

/  J 

(33) 


where 


22=4/«V 


(34) 


13 


and  n  is  any  positive,  non-zero  integer.  Neglecting  spatial  details,  it  is  interesting  to  note  the 
time  constant  r  of  the  exponential  tenn, 


r  = 


A  *  * 

1  *  2 

k  n 


(35) 


Next  we  considered  the  steady  state  solution  to  equation  32  with  boundary  conditions  of  zero 
temperature  flux  at  x*  =  0  and  T  =  0  at  x*  =  1 : 

c 


(36) 


Equation  36  represents  the  final  temperature  if  a  current  step  input  were  applied,  and  equation  35 
represents  how  quickly  that  temperature  is  attained.  A  hypothetical  step  response  for  x  =  0  is 


T * 


/ 2  / 

\( 

2k  V 


1-e 


-t  It 


(37) 


with  an  initial  slope  of 


2 k-  4 c>; 


7JV 

8c>: 


(38) 


For  short  times,  an  approximation  of  the  solution  is  a  linear  line  starting  at  x*  =  0  and  extending 
upward  with  a  slope  given  by  equation  38.  Thus,  we  expected  the  time  to  melt  to  be  affected  by 
changes  in  cv ,  pm  but  not  k ,  and  figure  8  demonstrates  this  expectation. 


We  note  that  the  use  of  time  for  the  central  node  to  reach  its  melting  temperature  may  not  be  the 
best  parameter  for  a  sensitivity  study.  Alternative  figures  of  merit  include  the  standard  deviation 
of  all  nodes  in  the  bridge  region  after  a  specified  time,  such  as  the  time  for  the  central  node  to 
reach  the  melting  point  or  the  time  for  the  average  of  all  bridge  region  nodes  to  reach  the  melting 
point.  Such  alternative  figures  of  merit  would  indicate  how  changing  input  parameters  might 
improve  or  impair  temperature  unifonnity,  which  is  expected  to  determine  post  melt  behavior. 

Of  course,  the  best  figure  of  merit  would  be  the  time  for  the  detachable  flyer  above  the  bridge 
wire  to  detach,  or  its  final  velocity  after  a  certain  distance,  but  such  calculations  are  the  subjects 
of  future  work. 
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4.  Conclusions  and  Future  Work 


We  have  modeled  the  pre-melt  behavior  of  an  exploding  foil  initiator  or  slapper.  The  results  in 
tenns  of  time  for  the  center  of  the  exploding  foil  to  reach  the  material  melting  temperature  was 
within  2%  of  a  similar,  published  modeling  study,  and  approximately  within  20%  of  the 
experimentally-observed  time -to-burst  for  specific  geometries  fabricated  at  ARL.  We  performed 
a  sensitivity  study  on  several  geometric  and  material  properties,  and  found  the  melting  time  to  be 
most  sensitive  to  bridge  wire  thickness,  and  all  material  properties  except  thermal  conductivity. 
We  used  simple,  analytical  arguments  to  explain  these  results. 

Of  course,  the  predicted  melting  time  of  one  location  in  the  bridge  wire  is  not  necessarily 
comparable  to  the  experimental  time  to  burst,  and  we  anticipate  two  future  technical  reports 
detailing  more  complex  models  to  treat  post-melt  behavior.  We  anticipate  the  next  report  to 
focus  on  the  post-melt  phase  of  the  metal  as  it  transitions  from  solid  to  liquid,  with  each  finite 
element  absorbing  the  latent  heat  of  fusion  and  with  material  properties  such  as  conductivity 
changing  according  to  empirically-determined  models.  The  liquid  regions  will  then  be  free  to 
expand  and  flow,  requiring  the  coupling  of  a  fluid  mechanics  code  with  the  existing  solid 
mechanics  model.  Portions  of  the  model  should  quickly  reach  evaporation  temperatures 
requiring  another  phase  change  into  a  gas.  Finally,  we  anticipate  the  last  technical  report  to 
cover  the  incorporation  of  a  surrounding  flyer  material  and  substrate,  which  will  allow 
accurately-modeled  heat  transfer  and  fluid-structure  interaction  as  the  bridge  wire  melts  and 
expands.  We  plan  to  model  the  non-linear,  plastic  deformation  of  the  flyer  material  either 
through  an  appropriate  code  allowing  a  fully-defonnable  mesh,  or  to  investigate  the  possibility  of 
modeling  the  problem  with  fluid  mechanics  using  a  hypothetical  fluid  with  a  large,  non- 
Newtonian  viscosity.  It  may  also  be  possible  to  model  the  structure  with  a  standard,  linear 
stress-strain  model,  and  to  simply  monitor  when  plastic  deformation  becomes  too  excessive  to 
trust  the  results. 

The  path  forward  will  result  in  a  predictive  model  allowing  any  number  of  sensitivity  studies. 
More  importantly  for  model  verification,  model  outputs  will  be  available  for  direct  comparison 
to  experimental  data  such  as  time  to  burst  and  flyer  velocity,  and  should  greatly  enhance  our 
ability  to  design  better  bridge  wire  devices. 
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Appendix.  AN SYS  Input  File  Listing 


/CWD, 'C:\Documents  and  Settings\christopher.morrisl7\My  Documents\projects\fuzing\ansys' 
!/CWD,'C:\Documents  and  Settings\christopher.morrisl7\My  Documents\Runspace’ 

/FILN  AME  ,fuze_2  d_45  tap_OpO  8  fdlet_Cu,  1 

!The  problem  is  as  follows: 
thk  =  2.5e-6  !m 
a  =  .le-3  !m 

10  =  3800*(a/thk)/15  !Amps 
t_I  =  370e-9  !s 
T  inf  =  300  !K 
eV_to_K  =  11604  !K/eV 

!  Non-dimensional  temperatures  at  which  to  evaluate  properties 
T1=0 
T2=0.05 
T3=0.1 
T4=0.5 
T5=l 
T6=5 


!  Thus,  the  non-dimensional  values  become: 
a_star=l 

I0_star=l  !  Current  input  amplitude 


brlen  =  astar 
brwid  =  astar 
in_wid  =  5*  as  tar 
outwid  =  inwid 
taper  =  45  !  degrees 
fillet  =  0.08*a_star 

brmat  =  1  !  l=Cu,  2=Au 

e_size  =  a_star/30 

ndivb  =  100  !  number  of  horizontal  element  divisions  in  the  bridge 
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BUILDING  MODEL 


/PREP7  [Enter  the  preprocessor 

!  First,  delete  all  previous  geometry,  if  it's  there... 

ALLSEL 

ACLEAR,ALL 

ADELE,ALL,„1 

ldele,all,„l 

kdele,all,„l 

!/NERR,0, 10000  ISupress  warning  messages 

/triad, off 

/PNUM,KP,0 

/PNUM,LINE,  1 

/PNUM ,  ARE  A,  0 

/PNUM,VOLU,0 

/PNUM, NODE, 0 

/PNUM,SVAL,0 

/PNUM,ELEM,0 

/NUM,0 


**  MATERIAL  PROPERTIES  ** 


TOFFST,0  !  Temperature  offset  from  absolute  zero  (don't  let  ANSYS  apply  any  other  offsets) 
EMUNIT,EPZRO,  1  !  I  don't  think  permittivity  is  used  in  these  calculations,  but  I'll  set 

lit  to  a  non-dimensional  value  of  1  anyway 
[(like  I  would  for  electrostatic  analyses). 

[Material  1:  Copper 

T  melt  =  1358  !K 

rho_m  =  8920  !kg/mA3 

!— Coeff.'s  from  Burgess  resistivity  model 

Cl  =  -4.12e-5 

C2  =  0.1 13 

C3=  1.145 

!  Coef  s  on  fits  for  Cv(T),  k(T) 

cv3  =  0 

cv2  =  0 

cvl  =  0.1004 

cv0  =  355.21 
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kl  =  -0.0685 
kO  =  420.75 


sigma_0=100000/(Cl+C2*((T_inf/eV_to_K)**C3)) 

UIMP,l,DENS,„rho_m/(  ((I0**2)*(t_I**3))  /  (sigma_0*(a**6))  ) 

Tdim  1  =T  1  *  (T_melt-T_inf)+T_inf 

Tdim2=T2  *  (T_melt-T_inf)+T_inf 

Tdim3=T3  *  (T_melt-T_inf)+T_inf 

Tdim4=T4  *  (T_melt-T_inf)+T_inf 

Tdim5=T5  *  (T_melt-T_inf)+T_inf 

Tdim6=T  6  *  (T_melt-T_inf)+T_inf 


!  Specific  heat 
MPTEMP,,,,,,,, 

MPTEMP,1,T1 
MPTEMP, 2,T2 
MPTEMP, 3, T3 
MPTEMP, 4, T4 
MPTEMP, 5, T5 
MPTEMP, 6, T6 
MPDELE,C,1 

MPDATA,C,l,,(cv3*Tdiml**3  +  cv2*Tdiml**2  +  cvl*Tdiml  +  cvO) 
/((a**2)/(t_I**2*(T_melt-T_inf))) 

MPDATA,C,l,,(cv3*Tdim2**3  +  cv2*Tdim2**2  +  cvl*Tdim2  +  cvO) 
/((a!|”|52)/(t_IHs*2*(T_melt-T_inf))) 

MPDATA,C,l,,(cv3*Tdim3**3  +  cv2*Tdim3**2  +  cvl*Tdim3  +  cvO) 
/((aH5H52)/(t_IH”|52*(T_melt-T_inf))) 

MPDATA,C,l,,(cv3*Tdim4*!|53  +  cv2*Tdim4**2  +  cvl*Tdim4  +  cvO) 
/((a**2)/(t_I**2*(T_melt-T_inf))) 

MPDATA,C,l,,(cv3*Tdim5**3  +  cv2*Tdim5**2  +  cvl*Tdim5  +  cvO) 
/((a!|”|52)/(t_IHs*2*(T_melt-T_inf))) 

MPDATA,C,l,,(cv3*Tdim6**3  +  cv2*Tdim6**2  +  cvl*Tdim6  +  cvO) 
/((a**2)/(t_I**2*(T_melt-T_inf))) 

!  Thennal  conductivity 
MPTEMP,,,,,,,, 

MPTEMP,  1,T1 
MPTEMP, 2,T2 
MPTEMP, 3, T3 
MPTEMP, 4, T4 
MPTEMP, 5, T5 
MPTEMP, 6, T6 
MPDELE,KXX,  1 
MPDELE,KYY,  1 

MPD  AT  A  ,KXX ,  1  „(kl  *(T  1  *(T_melt-T_inf)+T_inf)+kO)/((IO*  *2)/((T_melt- 
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T_inf)*sigma_0*a**2)) 

MPDATA,KXX,l,,(kl*(T2*(T_melt-T_inf)+T_inf)+k0)/((I0**2)/((T_melt- 

T_inf)*sigma_0*a**2)) 

MPD  ATA,KXX,  1  „(k  1  *  (T3  *(T_melt-T_inf)+T_inf)+kO)/((IO*  *2)/((T_melt- 
T_inf)*sigma_0*a**2)) 

MPD  AT  A,KXX,  1 ,  ,(k  1  *  (T4  *  (T_melt-T_inf)+T_inf)+kO)/((IO  *  *2)/((T_melt- 
T_inf)*sigma_0*a**2)) 

MPD  AT  A,KXX,  1  „(k  1  *  (T5  *  (T_melt-T_inf)+T_inf)+kO)/((IO  *  *2)/((T_melt- 
T_inf)*sigma_0*a**2)) 

MPD  AT  A,KXX,  1 ,  ,(k  1  *  (T6  *  (T_melt-T_inf)+T_inf)+kO)/((IO  *  *2)/((T_melt- 
Tinf)  *  s  igmaO  *  a  *  *  2  )  ) 

!  Resistivity 
MPTEMP,,,,,,,, 

MPTEMP,1,T1 
MPTEMP, 2, T2 
MPTEMP, 3, T3 
MPTEMP, 4, T4 
MPTEMP, 5, T5 
MPTEMP, 6, T6 
MPDELE,RSVX,  1 
MPDELE,RSVY,  1 

MPDATA,RSVX,  1  „  l/(  1 00000/(C  1+C2*(  ((Tl*(T_melt-T_inf)+T_inf)  /  eV_to_K)**C3) 
)/sigma_0) 

MPDATA,RSVX,  1  „  l/(  1 00000/(C  1+C2*(  ((T2*(T_melt-T_inf)+T_inf)  /  eV_to_K)**C3) 
)/sigma_0) 

MPDATA,RSVX,  1  „ l/(  1 00000/(C  1+C2*(  ((T3*(T_melt-T_inf)+T_inf)  /  eV_to_K)**C3) 
)/sigma_0) 

MPDATA,RSVX,  1  „  l/(  1 00000/(C  1+C2*(  ((T4*(T_melt-T_inf)+T_inf)  /  eV_to_K)**C3) 
)/sigma_0) 

MPDATA,RSVX,  1  „ l/(  1 00000/(C  1+C2*(  ((T5*(T_melt-T_inf)+T_inf)  /  eV_to_K)**C3) 
)/sigma_0) 

MPDATA,RSVX,  1  „ l/(  1 00000/(C  1+C2*(  ((T6*(T_melt-T_inf)+T_inf)  /  eV_to_K)**C3) 
)/sigma_0) 


■End  of  Material  props 


*AFUN,deg 

K,l,-br_len/2  -(out_wid/2-br_wid/2)/tan(taper),  +out_wid/2 

K,2,-br_len/2,br_wid/2 

K,3,br_len/2,br_wid/2 

K,4,br_len/2  +  (out_wid/2  -  br_wid/2)/tan(taper),  +out_wid/2 
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L,l,2 

L,2,3 

L,3,4 

LFILLT,  1,2, fillet 
LFILLT,2, 3,  fillet 

k,2,-br_len/2  -(out_wid/2-br_wid/2)/tan( taper), 0  !KP  #2  got  replaced  by  5,6 

*get,xcoord,kp,5,loc,x 

K,9,xcoord,0 

*get,xcoord,kp,6,loc,x 

K,10,xcoord,0 

*get,xcoord,kp,7,loc,x 

K,1  l,xcoord,0 

*get,xcoord,kp,8,loc,x 

K,12,xcoord,0 

K,3,br_len/2  +  (out_wid/2  -  br_wid/2)/tan(taper),0  !KP  #3  got  replaced  by  5,6 


L,l,2 

L,2,9 

L,9,10 

L,  1 0, 1 1 

L,  11,12 

L,12,3 

L,4,3 

L,5,9 

L,6,10 

L,7,ll 

L,8,12 


[Create  Areas 
AL,  1,6, 7, 13  !  1 

AL,8,14,4,13  !2 

AL, 9, 15, 2, 14  !3 

AL, 10, 16, 5, 15  !4 

AL, 11, 12, 3, 16  !5 


lsel,s,line„12,16 

lsel,a,line„6 
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[vertical  lines 


lesize,all,„ndivb/2, 1 0 
lsel,s,line„4,8,4 
lsel,a,line„5,10,5 
*if,nint(fillet*ndivb),eq,0,then 

lesize,all,„3  [fillet  lines 

*else 

lesize, all, ,,3*NINT(fillet*ndivb)  [fillet  lines 

*endif 

lsel,s,line„2,9,7 

lesize,all,„ndivb,-10  [horizontal  lines  in  bridge 
lsel,s,line„l,7,6 

lesize, all, ,,40,1/60  [horizontal  lines  leading  into  bridge 
lsel,s,line„3,l  1,8 

lesize, all,,, 40, 60  [horizontal  lines  leading  into  bridge 

allsel 


Select  areas  and  assign  attributes 


[Set  element  types 
ET,  1  ,PLANE223, 1 1 0 

[-select  membrane  areas 
ASEL,S, AREA,,  1,5 
AATT,br_mat,  1,1,0 

[  **  Mesh  ** 

ESHAPE,2  [Use  mapped  mesh 

ESIZE,e_size 
ASEL,S,  area,,  1,5 
AMESH,ALL 

allsel 


END  OF  BUILDING  MODEL 


t _ 

!—  SOLUTION 
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/solu 

ANTYPE, TRANS, NEW 
!  Apply  initial  conditions  to  all  nodes: 

IC, all, TEMP, 0 
IC, all, volt, 0 

KBC,0  !  ramped  boundary  conditions 
OUTRES, ERASE 

*afun,rad 

PI=ACOS(-l) 

tstep=  (l+T_inf/(T_melt-T_inf)) 

tstep=tstep  *  (c  v  1  *Tdim  1  +cvO)  *  (t_I*  *2)  *  (  T_melt-T_inf)/(a  *  *  2) 

tstep=tstep  *  rho_m  *  sigma_0  *(a**6)/(I0**2)/ (t_I  *  *  3  ) 

tstep=tstep/ 4/(PI  *  *  2) 

tstep=tstep**(l/3) 

tstep=tstep/30 

t  total =  .5 


t _ 

*do,ls,tstep,t_total,tstep 

time, Is 

/com, 

/com, 

/com, 

/com, 

/com,  incrementing  is 
/com, 

/com, 

/com, 

/com, 


Apply  BCs 


lsel,s, line,,  1,5 
nsll,s,l 
f, all, amps, 0 
f, all, heat, 0 
lsel,all 


nsel,s,loc,y,0,0 
f, all, amps, 0 
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f,all,heat,0  !  symmetry  BC 

*afun,deg 

nsel,s,loc,x,br_len/2  +  (out_wid/2 
d, all, volt, 0 
d, all, temp, 0 


br_wid/2)/tan(  taper) 


nsel,s,loc,x,-br_len/2  -(out_wid/2-br_wid/2)/tan( taper), -br_len/2  -(out_wid/2- 

br_wid/2)/tan(  taper) 

nsel,u,loc,y,out_wid/2,out_wid/2 

*get,my_b,node„mxloc,y 

*GET,  nodes,  NODE,  0,  count 

*afun,rad 

PI=ACOS(-l) 

Iin=T0_star  *  sin(2  *PI  *ls) 

f, all, amps, 0.5*Iin/nodes/(l-((out_wid/2  -  my_b)/(out_wid/2)))  !  Current  is  applied  so 

!  the  total  current  into  the  bridge  is  Iin.  The 
!  factor  of  0.5  accounts  for  symmetry.  The  other 
!  factors  account  for  not  applying  current  to 
!  "walls" 


f, all, heat, 0 


allsel 

TRNOPT,FULL 

LUMPM,0 

!* 

!OUTRES,ALL,l 
!  OUTRE  S  ,nsol,last 

solve 

/postl 

!SET,LAST  ! read  last  set 

nsel,s,loc,x,0,0 

nsel,r,loc,y,0,0 

*get,nodenum, node, ,num, max 
*  get, my  T emp  ,node  ,nodenum,temp 

*if,myTemp,gt,l,exit  !  Stop  after  the  melting  temperature  is  reached 

allsel 

/solu 

antype,,rest  !  restart  the  analysis 
*enddo 
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! 


-  END  OF  SOLUTION 


/postl 

!  interpolate  time  to  melt: 

SET,  NEAR„„LS-tstep 
*  get, T  last, node, nodenum,  temp 

tm_star  =  (l-T_last)*(tstep/(myTemp-T_last))  +  (LS-tstep) 
set, last 

/EXPAND,2,RECT,HALF„0.0000 1 

PLNSOL,  temp,,  0 

/AUTO,l 

/ZOOM,l,SCRN, 0.340024, 0.006081, 0.596081, -0.210017 

*if,l,eq,0,then 

PLNS,TEMP, 

ANTIME,50,0.2,  ,1,2, 0,0.5 
*endif 

/POST26 

!FILE,'fuze_2d’,'rst','.' 

/UI,COLL,l 
NUMVAR,200 
SOLU,  191  ,NCMIT 
STORE, MERGE 
FILLDATA,  1 9 1  „„  1 , 1 
REAL  VAR,  19 1,191 
!* 

nsel,s,loc,x,0,0 

nsel,r,loc,y,0,0 

*get, nodenum, node, ,num, max 

NSOL, 2, nodenum, TEMP,,  TEMP_2 

STORE, MERGE 

XVAR,1 

PLVAR,2, 

allsel 

FINISH 
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